Robustness of quantum memories based on Majorana zero modes 



(N 
O 

o 

Q 

0^ 



^ : 

Oh- 

> ', 

a : 

^3 ■ 



> 
OO 

(N 



X 



L. Mazza/'2'E1 M. Rizzi/ M. D. Lukin,^ and J. I. Cirac^ 

' Max-Planck-Institut fur Quantenoptik, Hans-Kopfermann-Strasse 1, D-85748 Garching, Germany 
^Scuola Normale Superiore, piazza dei Cavalieri 7, 1-56127, Pisa, Italy 
'Physics Department, Harvard University, Cambridge, Massachusetts 02138, USA 

We analyze the robustness of a quantum memory based on Majorana modes in a Kitaev chain. 
We identify the optimal recovery operation acting on the memory in the presence of perturbations 
and evaluate its fidelity in different scenarios. We show that for time-dependent Hamiltonian per- 
turbations that preserve the topological features, the memory is robust even if the perturbation 
contains frequencies that lie well above the gap. We identify the condition that is responsible for 
this feature. At the same time we find that the memory is unstable with respect to particle losses. 



Introduction. Quantum memories for the reliable 
storage of quantum states over long times constitute 
one of the the basic ingredients for most applications in 
quantum information science P, The so-called self- 
protected quantum memories (SPQM) 0,13 constitute 
appealing candidates for realization of quantum memo- 
ries. Here, the quantum state is stored in the low-energy 
subspace of a many-body quantum system, and the in- 
teractions among its elementary constituents provide the 
protection against perturbations. As a key feature, ac- 
tive error correction is not required. A single recovery 
operation at the time of retrieving the quantum state [5| 
is sufficient. The degree of protection may, however, not 
be as strong as that offered by fault-tolerant error cor- 
rection In recent years, a significant effort has been 
devoted to quantify the robustness of SPQMs against dif- 
ferent types of perturbations 

Topological superconductors (TSc) hosting spatially- 
localized zero-energy Majorana modes are attrac- 
tive candidates for a SPQM. Those systems exhibit a 
gapped quasi-degenerate ground state subspace whose in- 
ternal energy splitting closes exponentially with the sys- 
tem size of th e sy stem, even in presence of small local 
imperfections [l8|. Thus, the dephasing time of a qubit 
encoded in that subspace grows exponentially with the 
size of the system. This intriguing property is associated 
to the topological character of the TSc [l5|, |l8| , and thus 
it is referred to as topological protection of the informa- 
tion. Moreover, the fermion parity superselection rule can 
be used to protect the qubit against flip errors, as long as 
the environment does not introduce single particles in the 
system [l^. Although Majorana modes are hard to find 
in Nature, several proposals to engineer them have been 



put forward [201 121[ . The experimental observation and 



study of such modes are now actively being pursued [22 



The robustness of quantum memories based on Majo- 
rana modes has been analyzed under two assumptions 
(see, eg. [H, 1^ 2^ 24|): (i) the quantum state is stored 
in the quasi-degenerate ground state of the (eventually) 
perturbed Hamiltonian [25| ; (ii) a particular recovery op- 
eration is applied. Whereas those assumptions are justi- 
fied in some scenarios, this is not always the case. For 



instance, in cases when the qubit is not stored in the 
ground subspace of the Hamiltonian, either because the 
latter is not known, or it suddenly changes in time, all 
the eigenstates of the Hamiltonian are involved in the 
dynamics. Thus the robustness cannot only rely on its 
ground state properties [l^. Moreover, even if the stan- 
dard recovery operation fails to retrieve the original state, 
there may exist more general operations that could ac- 
complish this task. In this case the information about the 
stored state may still exist in the system, but cannot be 
extracted with the standard recovery. 

In this Letter we study the robustness of a quantum 
memory based on the Majorana modes of the Kitaev 
chain [17| to both Hamiltonian perturbations and inter- 
actions with the environment. Focusing on the scaling of 
the memory time with the system size, we first develop 
the formalism to identify the optimal recovery operation 
in the presence of perturbations and evaluate the opti- 
mal recovery fidelity. We next consider the storage of a 
qubit in the ground state of a given Hamiltonian, while 
the evolution takes place according to a different one. 
We find that, as long as the original and (instantaneous) 
perturbed Hamiltonians are within the topological phase, 
the memory time grows exponentially with the system 
size, even if the perturbation contains frequencies well 
above the gap. This implies that protection in this dy- 
namical scenario is not only related to the topological 
character of the ground state of the Hamiltonian, but to 
a dynamical condition which we identify. Furthermore, 
this provides us with an example of a memory that is not 
stable with the standard syndrome-based recovery oper- 
ation 12|, |27| , but can still be made robust by carefully 
choosing a more general operation. We also investigate 
the effects of temperature, the coupling to a reservoir in- 
volving particle exch ange, a nd the restriction to Gaussian 



recovery operations 

Setup. We consider a Kitaev chain of 2A^ Majorana 



modes, {cj}, with c| 
Hamiltonian is (l7| : 



and {cj,Ck}+ = 'i5j,k- The 
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At zero temperature, it has a topological phase in the 
regime J\ < 2. There, the ground state is quasi- 
degenerate due to the existence of (nearly) zero energy 
Majorana modes, 2 localized at the edges. They can 
be expressed as linear combinations of the {c}, and define 
a Dirac mode = 5 {rhi + 11112) ■ Due to superselection 
rules, we need two such modes to build a qubit, and thus 
we consider two Kitaev chains, with (Dirac) modes a and 
h. We construct a set of Pauli operators: 

(2) 

and define a qubit in the even parity sector: span{|0) = 
|vac),|l) ^atfetjvac)} [s^. 

The chains have non-trivial dynamics and (possibly) 
interact with the environment. Let us denote by TDt the 
decoherence channel that describes the time evolution in- 
cluding the action of the protecting Hamiltonian ([T]). It 
maps the state of the qubit onto the state of the two 
chains 31|. The chains are kept far apart, so that they 
do not interact with each other, nor with the same en- 
vironment. In such a case, we can simplify the analysis 
by assuming that the second Dirac mode, b, does not 
undergo any dynamics. 

Optimal recovery of the information. After the stor- 
age time i, we attempt to retrieve the qubit by apply- 
ing a recovery channel, TZt, to the chain and the extra 
fcrmionic mode, b. It maps back N + 1 fermionic modes 
onto the qubit defined in terms of ([2]). The composite 
channel Tt Tit ° T^t describes the complete action on 
the qubit. TZt has to be chosen such that 7* is as close as 
possible to the identity channel. This can be quantified 
with the recovery fidelity: 



FiUt)^ / d^l^{cp\Tt{\ip){ip\)\v), 



(3) 



where the integral is over the surface of the Bloch sphere. 
We find that the optimal fidelity is (see the appendix) : 



f: 



lOpt 



umxFiTZt) = 1 + l\\p+{t) - P-mtr; (4) 
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where p±{t) = Vt{pq,±), Pq,± = (1±(T^)/2, and || • ||tr de- 
notes the trace norm. If we restrict the analysis to Gaus- 
sian recovery operations. 



' Kt IS Gaussian o D 



r+(<)-r_(t)||op. 

(5) 

Here, T±{t) are the 2A^ x 2N covariance matrices (CM) 



itr[p±CjCfc] (see |28[), and 



of P±{t) defined as [Fj 

II • Hop is the operator norm. Note that Fq'I provides a 
lower bound to the optimal recovery fidelity. Expressions 
(HI and dl]) can be understood as follows. Due to the 
imperfections, the distance between any two orthogonal 
qubit states decreases with time, and no physical recov- 
ery action can increase it again [32| . Thus, the fidelity is 




.Fjj = 0.9985 



FIG. 1: Robustness of the quantum memory against Hamil- 
tonian perturbations: /xq = and: (a,c) /i_ = J, /i+ = 1.5J; 
(b) /i_ = 2. 5 J, 11+ = 3 J. Upper panels (a,b): Fidelity as a 
function of time for different chain lengtlis, N = 8, 16, 24, 32, 
40 and 48. Lower panel (c): Memory time to as function of 
TV for different values of the fidelity threshold Fq. We consid- 
ered = 101 realizations uniformly distributed in the range 
[/x_,/i+]. Data show convergence behavior for A^d — > 00 (see 
the appendix) and can therefore be considered as an approx- 
imation of the continuum situation. 



entirely characterized by the distance between the states 
evolved out of two orthogonal vectors [s^l. The recov- 
ery operation that reaches the optimal fidelity tries to 
align back the image of the Bloch sphere and consists 
of a unitary operation followed by the discarding of all 
the extra fermionic modes. In the general case, it may 
involve iV-body terms acting on all the fermionic modes 
of the chain and its specific form depends on the singular 
value decomposition of /5+ (t) — /5_ (t) (see the appendix 
for details). In the simpler Gaussian case, the unitary op- 
eration is an evolution with a Hamiltonian similar to the 
original one (see the appendix). 

Hamiltonian perturbations. We consider first the ef- 
fect of perturbations of the Hamiltonian: we assume the 
perfect encoding of the qubit using the two lowest-energy 
eigenstates of Hq ~ (/ig), and then we add an unknown 
perturbation, V. Such perturbation is randomly chosen 
from a set according to some measure Vy. We then add 



incoherently the corresponding evolutions 31 1 



Vtipq) = / di^y e-^(^"+^)*p(0)e»(^"+^)*. (6) 



We choose V quadratic in the Majorana operators, so 
that we can efficiently evolve each term of the integral 



34| . can be computed by using a discrete measure 



and taking into account that Vtipq) has support in a 
relatively small subspace (see the appendix). 

Let us start with the results for a perturbation such 
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FIG. 2: Robustness of the quantum memory against time- 
dependent Hamiltonian perturbations: /io = and /i_ = J, 
fi+ = 1.5J. The Hamiltonian time-evolution is swapped be- 
tween H{fj,-) and H{fj,+ ) every St — J~^/4. The figure shows 
the optimal fidelity as a function of time for different sizes of 
the chain, A^. We take = 101 realizations. 



that Hq + V = H{ij). Wc take a distribution of /i which 
consists of A^d discrete points equally spaced in /i+J. 
Figure [Ha,b) displays F°p* as a function of t for differ- 
ent system sizes. We have chosen ^± such that the per- 
turbed Hamiltonians lie inside (a) and outside (b) the 
topological phase, respectively. In the first case, the de- 
cay time of F°^'^ does not depend on the size N , whereas 
in the latter does. This is quantified in panel (c), which 
displays fqfjV), the time at which a prescribed fidelity 
threshold [l^, fo, is crossed versus N. We take a poly- 
nomial fit of F"^*" so that to(N) is then obtained at the 
interception of the fitting curve with Fq . Results are com- 
patible with an exponential growth of the memory time 
with A^. We have also considered inhomogeneous local 
perturbations obtaining the same qualitative results. 

To analyze time-dependent perturbations, we focused 
on a square wave, fj,{t) ~ p, + S/j, ■ sgn [sin(27ra;<)] where 
{p.,Sfi} = ± ^_)/2 is averaged over = 101 re- 
alizations. In Figure [2] we have taken u = 2 J, which is 
larger than the Hamiltonian gap (of the order of J) . Re- 
markably, even in this case the memory time grows with 
the system size N. 

These results indicate that the memory can withstand 
perturbations that keep the system within the topolog- 
ical phase, despite the initial state having plenty of ex- 
citations (in terms of each perturbed Hamiltonian). The 
stability can be traced back (see the appendix) to the 
fact that for any pair of Hamiltonians Ha and Hi, in the 
topological phase, 



{0\exp{iHat)exp{-iHbt)\Q) « {l\cxp{iHat) exp{-iHbt)\l) 

(7) 

where the difference decays exponentially with N (|0) 
and |1) are the two ground states for the original Hamil- 
tonian, Hq). This means that the excitations generated 
by any pair of evolutions have almost the same overlap 
(in modulus and phase) independent on whether they 
started, in |0) or |1). This appears to be a consequence of 
the topological nature of our system and therefore should 



FIG. 3: Optimal Gaussian fidelity as a function of time for 
different system sizes N: no = and /i_ = J, fi+ — 1.5J. 
We considered at most A^d = 101 realizations uniformly dis- 
tributed in the range [fJ,-, fJ.+]; plotted values for A'd — >■ cxj are 
obtained via scaling as 1/A^d (see the appendix). 
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FIG. 4: Upper bound to the optimal fidelity in presence of 
thermal modes, (left) Temperature is fixed: jS^^ — 2.5 J and 
the figure shows the upper bound as a function of time for 
different lengths of the chain TV = 8, 16, 24, 32, 40, 48. (right) 
Size is fixed: A'^ = 32 and the figure shows the upper bound 
as a function of time for different temperatures /3~^ = 0.5J, 
l.OJ, 2.5J, 5J and lOJ. We take A'^d = 101 realizations. 



likely be present in other systems as well. 

Time- independent perturbations have also been re- 
cently considered in [1^, obtaining different conclusions. 
In particular, in their analysis the memory time does not 
grow exponentially with the system size. The discrep- 
ancy originates from their choice of the syndrome based 
recovery operation (27| . which, according to our results, 
is not the optimal one. Nevertheless, they introduce ran- 
dom site-dependent variations of the chemical potential 
to localize the excitations and increase the memory time 
(see also [sf" 
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Our analysis shows that this is not 
necessary if the optimal recovery operation is used. 

We next examine memory performance for a restricted 
set of Gaussian recovery operations. Figure [3] displays 
the optimal Gaussian fidelity in the presence of the same 
perturbations as in Fig. [TJa): the fidelity improves with 
the size of the system but saturates for sufficiently large 
N. Remarkably, similar results are obtained when all the 
modes different from a (i.e. corresponding to non-zero 
Majorana modes) are initially in thermal equilibrium at 
a temperature as shown in Fig. ID Because we deal 
with convex combinations of mixed states, we can only 
plot an upper bound Fj^^ to the optimal fidelity (see the 
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FIG. 5: Effect-9 of particle losses according to equation ((Sjl; 
fio = and F = J. (left) Upper and lower bound to the opti- 
mal fidelity as a function of time, showing a clear exponential 
scaling, (right) Upper bound to the optimal fidelity as a func- 
tion of time for different system sizes: N — 8, 16, 24 and 32. 
Results coincide exactly and are indistinguishable. 



appendix). For temperatures above the gap, a clear sat- 
uration behavior is displayed, see Fig. |4] (left). For tem- 
peratures below the gap, our results are not conclusive 
(not shown) , but lowering the temperature while keeping 
a fixed value of N clearly increases the memory time, as 
shown in Fig. S] (right). Thus, the temperature defines an 
effective size up to which topological protection can oc- 
cur. These result regarding the Gaussian recovery opera- 
tion can be understood as follows: (i) although each term 
in the integral ^ is Gaussian, its sum is not, and thus 
the density operator is not Gaussian either; (ii) a Gaus- 
sian recovery operation can only depend on the covariant 
matrix, F, which for (pure) non-Gaussian states coincides 



with that of a mixed Gaussian state [38|. Therefore, re- 
stricting to Gaussian recoveries has the same effect as 
considering mixed states, and this explains the similarity 
between this case and that of finite temperature. 

Markovian interaction with a bath. We analyze now 
the interaction of the system with a bath with which it 
can interchange particles via a Lindblad master equation: 



N 



dtp = -i 



vY,[dnpdi~\{diL,p}+]. (8) 



Results in figure [5] demonstrate that the information is 
corrupted at a rate which does not depend on the system 
size. As expected, the exchange of single particles with a 
bath mixes states of different parity, and thus gives rise 
to irreversible losses. 

Conclusions. We have studied the possibility of stor- 
ing quantum states in a Kitaev chain in the presence of 
perturbations and decoherence. In particular, we have in- 
vestigated for how long the information is still present 
somewhere in the Hilbert space corresponding to the 
chain. We have derived a formula that quantifies the op- 
timal (memory) fidelity for both general and Gaussian 
recovery operations. For both, static and time-dependent 
Hamiltonian perturbations within the topological phase, 
we have gathered numerical evidence that the memory 
time grows exponentially with the system size. We have 
identified the condition ([7]) which is responsible for this 
result. For Gaussian recovery operations or at finite tem- 
perature, the memory time saturates with the system 
size. We have also confirmed that when particle losses 
are present, the memory time is independent from the 
system size. We believe that our results can be extended 
to other topological superconductors based on Majorana 



16| 



as 



zero energy modes, like the Px + ipy model [1 
well as that our methods could have applications in the 
study of fractional quantum Hall systems. In particular, 
it will be very interesting to investigate more in detail 
the general relation between the the memory protection 
against dynamical perturbations, as analyzed here and il- 
lustrated by condition [71) , with the topological character 
of the Hamiltonian responsible for the dynamics. 
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Here d„ = \ {c2n-i +ic2n) annihilates one fermion in the 
n-th physical site, and Hq is the protecting Hamiltonian. 
This equation transforms Gaussian states into Gaussian 
states, and thus can be rewritten in terms of the CM [s^ . 

Given previous works [l^, one would expect that the 
memory time does not scale with the system size. How- 
ever, since optimal recovery operations have not yet been 
considered in this context, we investigate whether the re- 
sult extends to the optimal recovery, and thus compute 
. Unfortunately, the trace- norm appearing in (j4|) can- 
not be directly computed in terms of the CM. In or- 
der to circumvent this problem, we bound F°^^ with 

hwp.At) - P.Amtr < {i- Fu(2x,+,px,-r)'/' i, 

where Fjj is the Uhlmann fidelity j6|, which for Gaus- 
sian states is a function of their CM (see the appendix) . 



Appendix: Fermionic Gaussian States 

The formalism of fermionic gaussian states (FGS) [2^ 
is particularly useful in the treatment of quadratic 
fermionic theories, which include a wide class of topo- 
logical materials Not only FGS comprise the 
ground states and thermal states of such Hamiltoni- 
ans (via Bogoliubov transform) but also they describe 
dynamical evolution under some master equations (ex- 
actly) [1^ or in presence of moderate interactions (ap- 
proximately) [4^. For a N modes system, calculations are 
restricted on a space scaling only as N rather than expo- 
nentially. Indeed, these states are fully characterized by 
the sole (antisymmetric and real valued) covariance ma- 
trix (CM), i.e. the expectation values of quadratic com- 
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binations of fields, whereas all higher moments can be 
deduced via Wick theorem 41 1 . 

Given N fermionic modes, we can conveniently rewrite 
the 2N canonical Dirac fermionic operators {d|^^}i=i...Ar 
(with {di, dj} = 0; {di, d^j} = Sij) in terms of Majorana 
operators, i.e. fermionic operators which are real, Hermi- 
tian and unitary; 



Cj,2 = --i{dj - rf]); {cm, c„} = 2(5„i,„, 

(9) 

where m e {(1, 1), . . . (A^, 2)} glues two sub-indices to- 
gether for brevity. By a bit of simple algebra, Eq.(l) of 
the Letter is seen to be a rewriting of 

N-l ^ N 

h.c.) -/x^djd;.. 



AdU^ 



where J = A = 1. Canonical transformations can be 
represented by orthogonal real matrices, Ck — >■ c'j, — 
^[Ok.ici, as well as by a unitary rotation c'j^ — UckU'^ 
in Fock space in the case of O € 0{2m), i.e. 

detO = 1, the relation reads O = exp(A) ^ U = 
exp (—Aa^pCaCp/A) apart from an arbitrary phase. The 

number parity operator reads P = (— = 
Ylf^ Ck and is almost invariant under canonical trans- 
formations: P' = dot O ■ P. 

A TV-modes FGS is a A^-fermions state which has a 
density operator of the form p = 11^=1 Pa, with 



1 + e-'B^ 



,l) 



1 — iXaCa iCa 2 



2cosh(/3Q/2) 



(10) 

where Aq = tanh(/3Q,/2), and the Oq ■* (cq^o-) are the eigen- 
modes of the density operator. One can easily verify that 
Trp = 1, whereas Trp^ = + Xl)/2, i.e. the state 

is pure if and only if = ±1; moreover, p is positive ~ 
and thus a well-defined density operator - if and only 
if Aq G [—1,1]. FGS automatically satisfy the super- 
selection rule of the fermionic parity P, and therefore 
their density matrices can be expressed as a direct sum 

P — Pcvcn ffi Podd- 

The skew-symmetric covariance matrix (CM) of a FGS 
p is defined as the table of quadratic expectation values: 



- m , n — 2 Tr [p ( Cm Cn 



0]: 



e 



-Xa 




(11) 

where the second expression is given in the eigenbasis of 
Ea. ([TU|) . Under a canonical transformation O, the CM 
transforms as F' — OTO^ . The CM completely charac- 
terizes the properties of a FGS, as elegantly stated by the 
following reformulation of the Wick's theorem [4l|: 



Tr[p Co 



Pf 



Fl 



lQl...Q2p 



(12) 



where Fl is the restriction of F to the modes 

lQi...Q!2p 

{ai . . .a2p}, and Pf denotes its Pfaffian. This allows to 
simulate FGS efficiently with classical computers. 



The squared overlap of two FGS p and a is [28| : 



Tt[p&] 




V F 



(13) 



Moreover, generalizations of formula [T3] can be derived 
for any two Gaussian operators, as explained in [28|. Be- 
cause the time-evolution U (t) under a quadratic Hamil- 
tonian _ff is a Gaussian operator, such formulas for 
TT[pU{t)] and TT[U'{t)U"{t)] have been widely used in 
the Letter. Moreover, the CM of U (t) p{0)U (t)^ is: F(t) = 
O{t)r{0)O(t)^ , where 0{t) = e"^* and T is the real 
skew-symmetric matrix such that H ~ ^ ^ Tj ^CjCk- 

Also the Uhlmann fidelity among mixed Gaussian 
states, Fu{p,a) = 'Yy\J p^l'^ap^l'^ , can be efliciently com- 
puted via their CM's. By using Eq. (fTO|) . we can define Hp 

such that /5^/^ = exp(— iJp)/y^FYexp(— 2^p) and the cor- 
responding imaginary-time evolution of the state a [i^ : 

p^''^ (jp"'"'''^ = Tr[/5(T] ■ (T7(r = 1); oiij^ — 



Tr 



(14) 

Since (T/(t) is still a Gaussian state, whose CM can be 
efficiently computed [iOl, the trace of its square root in 
Fu can be calculated by looking again at Eq. (fTO|) as 
above. 



Appendix: Optimal Recovery Operation 

We present the derivation of formula (4) of the Letter, 
which describes the fidelity of the optimal recovery op- 
eration acting after a given decoherence map, rewritten 
here for reading convenience: 

^ ■ ^1Px,+(t)-A.,-Wlltr; (15) 



-lopt 



( — maxF(7^t) 



6' 



where Px,±{t) = A(*x,±) and *:r,± = (1 ± CTi)/2. 

We first show that the fidelity of any recovery oper- 
ation is upper bounded: F{TZt) < WTZt and then 
construct an explicit recovery operation, 7?,°''*, which 
achieves the upper bound. 

Before presenting the proof we recall that for any 
bounded Hermitian operator X the trace norm is 
llXjltr = max^ tT{HX), and the maximization is re- 
stricted to Hermitian operators fulfilling ||if||op < 1, 
where ||>?||op denotes the maximum singular value of X. 

Moreover, F{TZt) is an average over the Bloch surface 
(see Eq. (3) in the Letter) which is best expressed as [33| : 



(16) 



We derive the upper bound. Because of the contractiv- 
ity of the trace distance under TZt 0i the inequality 

trKrtia'J] < ll'^«l|op||7;(a:,)|ltr < 1 • ||A(^L)lltr 



6 



holds, from which we obtain: 

FiTZt)<l + ^ J2 \\Po.At) - Pc.-mtr- (17) 



12 

Let us now specify (fT7|) to the case of Vt not acting 
on the fermionic mode Aq. The operator = aa^ — a^a 
is Hermitian and ||-ffz||op < 1; moreover tr[Hz Pz.±{t)] = 
±1. Thus we have: 

2 = tr[ff,a-,+ (i)-/5..^W)] < 

< WpzAt) - Pz,-{t)\\tr < 11*.,+ - *.,-||tr = 2 



SO that \\pz^+{t) ~ pz-{t)\\ti = 2. Furthermore, since 
Vf does not act on the anciha, we can write Py^rit) 
%yPxAt)y^y^ where %y = g-^^'^/^. Therefore: 



- p.,-(t)||tr = - Py,-{t)hr (18) 

The bound F(7^^) < F°^^ \/TZt foUows from the combi- 
nation of (|17p with these considerations. 

The recovery map 7?,°^' which achieves F{TZt) — 
is in the form: 



and the operators Ha are such that ti[Ha'Dt{aa)] = 
||2?t((T^)||tr- Hz has aheady been defined. The Ha can 
be interpreted as the observables to be measured in p{t) 
in order to reconstruct pq 
Let us rewrite 

p^,+(i) - p^,-{t) :^aR + R^a^; R = Pt(-5) 

We write the most general Hermitian operator: 

Hx = aSi + S'Jflt + a^aS2 + aa^ 

which must satisfy H^Hx < 1- We get 



tr 



Hx(aR + i?^at 



= tr 



SiR^ + RSl) at 



Using the left polar decomposition R = PU, where P ~ 

\/Mt is positive semi-definite and U is unitary, we have 
that the maximum is attained when Si = U, and S2 = 
S3 = 0. Therefore, the maximum is achieved by: 



Hx 



Furthermore, since both a and R change the fermionic 
parity, Hx,y do not. They also fulfill: 



tr 



HaPa.±{t) = ±^\\Px,+{t) - Px.~{t)\\ti-; a^x,y 



from which tr[HaPt(o-^)] = ||Pt((T^)||tr follows. 



The optimal recovery map (|19p is linear, trace preserv- 
ing, and it also preserves the fermionic parity, since Ha 
do. For it to be a valid quantum channel we have only 
to show that it is completely positive. We construct a 
unitary operator acting on all the fermionic modes, W, 



such that 7?.°''* (p) = tr 



, where the trace is taken 



WpW^f 

over the fermionic degrees of freedom which are not part 
of the qubit. The operator is: 



W 



1 



E 



^Hn 



(20) 



Using that H^ ~ 1 and HaHp 



Hy, where ea. 



— is the Levi-Civita symbol, one can show that W is unitary 



and that it satisfies 7?.°^' (p) = tr 



F{TZ°^^) saturates the bound in ([4]). 



WpW"^ 



Furthermore, 



Appendix: Gaussian Recovery Operation 

We present the derivation of formula (5) of the Letter, 
which describes the fidelity of the optimal Gaussian re- 
covery operation acting after a given decoherence map, 
rewritten here for reading convenience: 



TTiOpt 

= max 

' T^tis Gaussian 



- x,+ 



it) 



-(i)llop 



The proof is similar to that in the previous section. 
In this case, we need to explicitly consider the fact 
that the qubit is composed of fermions. In particular, 
we have to consider what happens if the decoherence 
channel changes the parity of the state. This occurs, 
for instance, when decoherence is caused by the inter- 
change of particles with a reservoir, as described by 
^t(p) = {b + b^)p{b + b^)- The final state contains all 
the information about the initial state of the qubit, even 
if now it has odd parity. The problem, that could be 
overlooked in the previous section, arises in this section 
because the Gaussian recovery operation cannot act in- 
dependently on the different parity sectors of the state, 
whereas a general recovery operation can. 

We denote rhi , mi the zero-energy modes of one Kitaev 
chain and 7713 and TO4 the modes of a decoherence-free 
ancilla. We express the Pauli operators of the qubit as: 



{a^ - a)(&t + &) = ic3C2; 
—i[a^ + a) (6'^ + 5) = icaCi; 



aa' 



(21) 
(22) 
(23) 



where the modes b^^'^ and o^'^^ have been defined in the 
Letter. Note the difference with the ct^ in Eq. (2) of the 
Letter. Since the a'a act on both the parity sectors of the 
qubit, they allow the recovery operation independently 
of whether the parity has been changed or not. 
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We recall that the action of a general Gaussian channel 
transforms a Gaussian M-modes state with CM F into a 
TV-modes state with CM V = BTB^ + A (see B 
and A are 2N x 2M and 2N x 2N matrices chosen such 
that Q, defined as: 



Q = 



A B 

-B^ 



(24) 



satisfies Q^Q < 1; ^ must be skew-symmetric. 

We denote TZct the Gaussian recovery operation, and 
Tc^t = T^G,t ° T^t- Moreover, we define: 



^out 



^rG,t(*o,+) ~ ^rG.t(*„.-) 



(25) 

^BtzA^bI, (26) 



the difference of covariance matrices corresponding to the 
states after the decoherence channel and after the recov- 
ery operation, respectively. Note that the matrices Aq, are 
2N X 2N matrices, whereas the matrices A°"* are 4x4 
matrices. The assumption that the decoherence channel 
does not act on the first two Majorana modes ci^2 is re- 
flected by some properties of Aq, which are best expressed 
considering the block structure: 



A. 



K -LP 



(27) 



where K', L and K" are 2x2, 2{N - 1) x 2, and 2{N - 
1) X 2{N — 1) matrices, respectively. We obtain: 

K'z = ^2 J; Lj; = ^/j, ; Ly = L^J; ^ (^i q 

(28) 

where li^2 are column vectors. Additionally, Lz, K'^ y and 
K'^ y are zero matrices. Thus: 



Aa;||opt — ||A„||opt — ll-^xllopt 5: 2 



(29) 



We can now show that -F(7^G,^) < F°^^ V7^G,t. The 
starting point is equation (|16p . modified as follows: 

^G,* = ^ + i^ E ^^KrGA<^'c.)]- (30) 



a—x,y,z 



Noting that the ct" are quadratic in the Majorana oper- 
ators and recalling the definition of CM, we obtain: 

tr KTgA^L)] = (Ar )^,,,, < II Ar llopt, (31) 

where (/3i,;02) = (3,2), (3,1), and (1,2) for a = x,y,z, 
respectively. The most general Gaussian recovery opera- 
tor 7^G,i yields: 



lArilop 



|B7jA„B^|lop < ||A„| 



(32) 



given that Q'^Q < 1. Since 2 > ||A.||op > (A^)i,2 = 
iKz)i,2 = 2, we get that HA^Hop = 2. Using we 
obtain the desired bound. 



We now provide an explicit Gaussian recovery oper- 
ation which attains the bound. The recovery operation 
consists of the application of a Gaussian unitary opera- 
tion Wq to the the system and in subsequently tracing 
out N — 2 modes of the system. To define Wq , consider 
the singular value decomposition of =■ UYA/"^ , where 
\J {V) is a unitary 2(A''-l)x2(7V-l) (2x2) matrix and 
S is a 2{]>S — 1) X 2 matrix. Clearly, it is also possible to 
construct U' and V such that = U'YJV and E' has at 
most two elements different from zero, Si_2 > ^2,1 > 0, 
named the singular values of L^. We define Wq to be 
the unitary transformation which is represented by an 
orthogonal transformation V' : 



Wg Iwl. = {V © U') 



(33) 



Physically, Wq, rotates all the information between an- 
cilla and system into the first two modes of the system. 
The other modes can be now traced out. Summarizing: 



7^°P',(/5) = tr WqpWI 



(34) 



The CM of Tl°^\{p) is a function of T, the CM of p: 



V*(p) ^ [(^'® f^')r(r^©c/'^)]|(,_4),(,_ 



-4) 



(35) 



Finally, let us prove that F{JR^q^^) saturates the boimd. 
Denote = T^q \ ° T^t ■ Clearly, 



'3.2 



-(Ar)3,,+(A^ 



a—x,y^z 



By construction of Wq, {A"^*-)^^ — l|Aa;||op. Since 
J commutes with every 2x2 orthogonal matrix, 
(A-')3^^ = ||A,||op = ||A,|iop. Finally, (Ar')i,2 = 
{Az)i2 — 2 because the orthogonal transformation V' 
leaves the covariance matrix of the fermionic mode Aq 
unchanged. Together with equation ([30]) . this shows that 
the recovery operation 72.^'^^ saturates the bound in ([5]). 



Appendix: Some Details on Hamiltonian 
Perturbations 

We study the effect of the disturbance described in 
Eq. (6) of the Letter; we use formula For simplicity, 
let us discrctize the integral (Aj/)^^/di^ ^cT^ Sf'' 
and define the x TVj matrices Gq and Gi: 



[Gr]j.k = (.9t| 



i H -i t ^ — i Hi- 1 I 



H,=Ho + V,; (36) 



with notation taken from the Letter. We now show that: 

i \\p,,+ {t) - /5.,-(i)llt. = {VG^d, VGM) (37) 

where (•, •) is the Hilbert-Schmidt inner product for ma- 
trices: {A, B) = tr [i?^^] . It follows that perfect recovcr- 
ability, namely ||px.+ (t) — px^-{t)\\^^ = 2, is equivalent to 
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FIG. 6: Dependence of F°^^ on A^d- These data refer to the 
topological perturbation of figure 1 of the Letter; A*' = 32, 
fj,o = 0, fi- — J, fi+ — 1.5J. The scaling is shown for five 
times. From up to down: 24J"\ 54J"\ 174J~\ 354J~\ 
474J~^. Data show a convergence behaviour for A''d — > oo. 



FIG. 7: Dependence of Fq^*" on These data refer to 

the topological perturbation of figure 1 of the Letter: A*" = 
32, /io = 0, ^- = J, /i+ = 1.5J. The scaling is shown for 
five times. From up to down: 12J"\ 24J"\ 54J"\ 174J~\ 
354J~^, 474J"^ Data show a A''^^ dependence: the A''d -t- cxd 
value is extrapolated with a linear fit (thin red lines). 



Go = Gi, i.e. to the two topological sectors of the theory 
"behaving" in the same way. 

First, a simple inspection shows that Ga = Gj., Go- > 
and II \/G^J7^!Ihs = 1- Let us consider the states |±) = 
{\go) ± \gi)) /V2 and use them to define the overlap ma- 
trix M : 



M 



with Af 



(r,r') 



A/(+-+) M(+^-) \ 1 / Go + Gi Go - Gi 
M(-,+) m(-.-) ) ^ 2\Ga-Gi Go + Gi 

(38) 

±j — \/|c. - o 'I''"')- The second equal- 
ity follows from the assumption of a decoherence-free 
mode. Let us now consider the set of 2A'd states: C = 
{e-'"^*\+)}f:^^ U {e-'"^^-)}f^i and an orthonormal 
basis B = such that: span S = span C A 

matrix Y representing the basis change: e~'^''*|+) = 
E,n:>,) and e-^^^'H = E^Y^+kJ^'^) defined 
by yyl' = M and can be computed via the eigenvalue 
decomposition of M. Given V unitary and D diagonal 
matrix such that: M = VDV\ a Y defined as y = V^/D 
is one such possible basis change. A simple algebra leads 
to: 



I 

-I 



y 



(39) 

We can obtain a more explicit expression of Y considering 
the eigenvalue decomposition of Go-; Ga = VaDaVj and 
observing that the matrix; 



V 



Vq -Vi 



(40) 



diagonalizes M. Some simple algebra shows that the sin- 



gular values of y^ 



I 
-I 



y are two-fold degenerate 



and coincide with the singular values of vGqGi . Equa- 
tion ([37|) follows from the positivity of Go- 

The numerical computation of Go- is efficient when Hj 
is a quadratic Hamiltonian so that Uj^k =i= e^^^'e"*^''* 
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FIG. 8: Singular values of the matrix in equation (|4ip . (Left) 
The quench Ifamiltonians are in the topological phase: /io = 0, 
fij — J, pk = 1.5J. (Right) The quench Hamiltonians are not 
in the topological phase: no = 0, fXj — 2.5J and /is = 3.0J. 



is a Gaussian operator (po = l.'7o-)(5o-| is FGS) and thus 
[Gcr]j,k is a function of the CM of and Uj^k- We warn 
the reader interested in reproducing the data that some 
care is required in order to obtain the proper phase. 

We checked the possibility of using the data to obtain 
information regarding the limit Nd — >■ oo. In figure [6] 
we show the dependence of F°p* for different values of 
Nd] data show a clear convergence behaviour, even if the 
functional form of such scaling was not found. In figure [7] 
we show the dependence of on Nd ; data show a clear 
A'd scaling, which allowed us to take the limit Nd oo 
with a linear fit. 

In the Letter we provide numerical evidence that the 
memory time of the system increases exponentially while 
letting the system size N ^ oo, which is equivalent to 

AT = Go - Gi ^^°°> 0. We first observe that: 



|Aj-fe| oc 



pf 



(41) 



where To- is the CM of p^ and T is the CM of Uj^, we 
assume tr[J7j.fe(t)] ^ in order to avoid lengthy regu- 
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larized expressions. The proportionality factor between 
l.h.s. and r.h.s. is bounded by 1. Direct numerical inspec- 

because 



7V- 



-> e 



-N 



tion of the matrices shows that Xj^k 
^(Fq + Fi) + T /tr[Uj^k{t)] has two singular values which 
scale exponentially to zero (see figure [5]). This is clear 
when Uj^k = I and T = since Fq + Fi has two zero 
eigenvalues corresponding to the zero-energy modes of 
the Hamiltonian Hq. When Uj,^ 7^ I it would be tempt- 
ing to interpret T/tr[C/j ^(i)] perturbation and in- 
voke some topological stability argument; unfortunately 
T/tr[C/j_fc(<)] is neither bounded nor it is short-range. 
Moreover, one would like to have an explanation which 
distinguishes the situation in which Uj.k{t) is the prod- 
uct of two time-evolution according to topological and 
non-topological Hamiltonians. Intuitively, the argument 
must reside on the fact that any topological Hamilto- 
nian spreads the localized zero-energy modes exponen- 
tially slower than non-topological Hamiltonians do. 

In the Letter we present results also for thermal states. 
Because in this case the matrix p±{t) is not a con- 
vex combination of a limited number of known pure 
states, we cannot exactly compute F°^^ . Using the fact 



that px -i 



1 Y^JVd 
ATh ^.7 = 1 



iH-it ; 



/5z:,±(0)e J and denoting 



p'^J^j. = e '^^*p2:,±(0)e'^^*, we compute the following up- 
per bound: 



IAt,+ - AT,-I|tr < min 

Treojv, iVd — II 



0) _ c,{^U)) 

rx,— 



< 



3 



The minimization over the set of permutation of Nd el- 
ements can be done with the so-called "Hungarian al- 
gorithm" or "Kuhn-Munkres algorithm" 4^- The com- 
putation of the Uhlmann fidelity Fu between two mixed 
FGS has been explained in the section on FGS. 

The fidelity of the optimal Gaussian operation requires 
the computation of the CM of Px,±{t), which are not 
FGS. Because these states are convex combination of 
known FGS and because the CM is a linear function of 
the state, the CM of Px.±{t) is: 



- x,± 



U) 

x,± 



(42) 



where v'^J^^ is the CM of p^;!^^. 
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